Thermodynamics and magnetic field profiles in \ow-k type-II superconductors 
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Two-dimensional low-K type-II superconductors are studied numerically within the Eilenberger 
equations of superconductivity. Depending on the Ginzburg-Landau parameter « = A/£ vortex- 
vortex interaction can be attractive or purely repulsive. The sign of interaction is manifested 
as a first (second) order phase transition from Meissner to the mixed state. Temperature and 
field dependence of the magnetic field distribution in low-K type-II superconductors with attractive 
intervortex interaction is calculated. Theoretical results are compared to the experiment. 
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Superconductors are divided into two groups, type-I 
and type-II, based on their behavior in external mag- 
netic field. The underlying physics is especially trans- 
parent within the Ginzburg-Landau (GL) model. It can 
be shown that the difference in thermodynamics origi- 
nates from the sign of vortex-vortex interaction, which 
can be attractive or repulsive. The sign of the inter- 
action is controlled by GL parameter k = A/£, ratio 
of penetration depth and coherence length near T c . In 
type-I superconductors (n < l/v2) vortices attract each 
other. In type-II superconductors (k > 1/ , vortex- 
vortex interaction is repulsive and the phase transition 
from Meissner to the mixed state at H cl is of the second 
order. 

It has been verified experimentally since a long time 
that away from T c even in type-II superconductors with 
k ~ 1 interaction between vortices can be attractive. 
Vortex attraction is manifested as a first order phase 
transition from Meissner to the mixed state. This type 
of superconductivity we denote as type-IIa to distinguish 
it from the type-lib superconductivity with purely re- 
pulsive vortex-vortex interaction. Obviously attractive 
interaction in type-IIa superconductors will turn into re- 
pulsive below some critical intervortex distance thus pre- 
venting vortices to collapse into the normal state as in 
type-I superconductors. Evidence for jump in magneti- 
zation and mixture of vortex and Meissner domains near 
H c i is found from many magnetic, calorimetric, decora- 
tion and neutron scattering measurements (experimental 
results are partially reviewed in Ref.EJ). Here we refer jusjk 
to a nicely performed experiment by Auer and Ullmaiercl 
in which all three types (type-I, Ha and lib) of supercon- 
ductivity were reproduced in TaN by tuning the value of 

K. 

Although low-K type-II superconductors are of less in- 
terest from the point of view of potential application, 
attractive vortex interaction is still a challenging prob- 
lem not fully understood. Most of the theoretical ef- 
forts to understand interaction of vortices as a function 
of parameter k relied on the extended GL model. Uti- 
lization of Bogomolnyi method to treat GL equations for 
k close to l/VSmgave the further impetus for this re- 
search direction™! However, range of applicability of GL 
formalism is confined to the region near T c . Therefore 



to calculate the full temperature dependence of k c (T), 
boundary between type-IIa and type-lib superconduc- 
tivity, one should rather turn to the microscopic theory 
of superconductivity. An approximate curve k c (T) has 
been calculated by KramerQ based on the asymptotic so- 
lution of Eilenberger equations of superconductivity. It 
is argued that the calculated curve is just the lower limit 
for the true one. BrandtQ has obtained jumps at H c i by 
variational approximation of his free energy functional. 
The only attempt to calculate numerically k c (T) with- 
out approximations besides those immanent to numerical 
procedure is due to Kleincx In this thorough study clean 
and dirty isotropic 3D superconductors are analyzed. 

The purpose of this work is two fold. First, we applied 
a new method to solve Eilenberger equations for a vor- 
tex lattice. It appears to be efficient and faster than the 
so-called "explosion method"ErQ. Second, we extended 
previous numerical study of 2D superconductorsa to the 
low-n case where the interesting attractive vortex-vortex 
interaction occurs. Though the quantitative results will 
be different from the 3D case, the qualitative (k, T) phase 
diagram is quite similar. This is explicitly shown by 
calculating phase diagram (k, T) of type-I, type-IIa and 
type-lib superconductivity. We focused on the magnetic 
and temperature dependence of calculated magnetic field 
profiles and compare the result with experimental data. 

Eilenberger version of BCS theory of superconductivity 
is proven to be a good starting point for both analytical 
and numerical calculations. For the isotropic gap Eilen- 
berger equations are 



[u, + u(V + iA)]/ = %, 



[w-u(V-iA)]/t = ** fl . 
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These are supplemented by the self-consistency equations 
for the gap function \& and vector-potential A 



*lnt = 2t^ 



w>0 



V x V x A 



2t 



Im ( u 9) 



(3) 
(4) 



2 



Equations are written in following dimensionless units: 
order parameter is measured in units ttT c , length in 
units R = Uo/(2ttT c ), vq is Fermi velocity, magnetic 
field in units Hq = $o/2irR 2 , vector potential in units 
A a = <P /2nR, energy in units E = (t:T c ) 2 N(0)R 3 . 
Eilcnberger parameter k is 

It is related to GL parameter n via k 2 — (7C(3)/18)k 2 in 
3D case and 

<«> 

in 2D case. Here ( is Riemmann's zeta function. Here, 
lu = t(2n + 1) is Matsubara frequency with integer n, t = 
T/T c is reduced temperature, u is unit vector directed 
along Fermi velocity. Eilenberger Green's functions /, 
P and g are normalized so that g = y/l — //t. Aver- 
age over the isotropic cylindrical Fermi surface < . . . > 
reduces to (1/2tt) J ...dip, average over polar angle p>. 
Expression for the free energy density difference between 
superconducting and normal state is given by, simplified 
with the help of Eilenberger and self-consistency equa- 
tions 



F = « 2 (Vx A)2-<]T/— £(*•/ + */t)\ ( 7 ) 

We use the following notation for spatial average C = 
(B/2ir) J ce nCdS, where B is magnetic induction. Apart 
from the extreme cases H ~ H c \ and H ~ i/ C 2, which 
allow for an analytical solution, the only tool at hand is 
numerical analysis. Standard way to solve Eilenberger 
equations for a vortex lattice is so-called "explosion" 
method (for details see RefsMD). Here, another approach 
which takes the advantage of periodicity of the vortex 
lattice has-. been adopted. The same approach is used 
by Brandto to solve numerically GL equations. For this 
purpose it will be useful to introduce auxiliariz-functions 
a and b through the following transfornmtioriiO 

2acxp(— i 2&exp(i0) 1 — ab 

1 + ab ' 1 + ab ' 9 ~ 1 + ab' 

where <f> is the phase of the order parameter, \& = 
\^\ exp (—i<fi). Equations for auxiliary functions a and b 
are decoupled and have the form of Riccati's differential 
equation 

uVa = -(w + iuA s ) a + ^ (l - a 2 ) , (8) 



uV6= (w + iuA s ) 6-^(1 -b 2 ). (9) 



In new gauge vector-potential A s = A — V<^> is propor- 
tional to the superfluid velocity. It diverges as 1/r at the 
vortex center (we put index s to denote its singular na- 
ture). Functions a and b are not independent but there 
is a relationship between them 

a*(r, — u, u) = b(r, u,u)) . (10) 

Therefore it is sufficient to solve only the equation for 
auxiliary function a. Also, we don't need to solve equa- 
tion (||) for all Fermi velocity directions because 

a(— r, — u, u>) = a(r, u, u>). (11) 

Reflection at the y axis will change the direction of the 
magnetic field. When combined with time reversal it 
leaves Hamiltonian of the superconducting system invari- 
ant. One can show, at least numerically, that the follow- 
ing equation holds 

a(x,y,-u x ,u y ) = a*(-y,-x,u x ,u y ). (12) 

Therefore it is enough to solve Eq. (|J) for velocity di- 
rections < ip < n/2 to recover a and b for all u via 
symmetry relations (|l0l), (|TT1), and (|l2"|). 

Singular vector-potential has the following property 

V x A s =B(r)-2ire z J2t(r-r i ) ) (13) 

where r, are positions of vortices and 6(r) is 2D Dirac 
function. For an ideal vortex lattice z-component of the 
magnetic field B(r) = V x A and moduo of the gap 
function |^| has the periodicity of the vortex lattice. The 
same is not true for ^ and A. However, superfluid ve- 
locity, i.e. A s , is invariant under the translation from 
one to another vortex lattice cell. Therefore one must 
assume that auxiliary functions a and b are also periodic 
functions. Then we expand a in Fourier series 

a(r) = ^a(Q)exp(*Qr), (14) 
Q 

where summation goes over reciprocal vortex lattice vec- 
tors Q. Fourier transform of a periodic function f(r) is 
defined as 

if(r)] FT = HQ) = J f(r) exp (-iQr) d 2 r. (15) 

cell 

We split the local magnetic field B(r) into constant part, 
magnetic induction B, and periodic part B'(r) with zero 
average over the vortex lattice cell. Similarly, vector- 
potential is splitted into two parts 

B(r) = B + B'(r) = V x A + V x A'. (16) 

We will assume that A'(r) is periodic with V • A' = 0, so 
that all uncertainty of the gauge is in A. Therefore.A' 
is periodic solution of the self-consistency equation (||) 

Of 

V 2 A' = 4lm£(u ff ). (17) 

rC 
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From the equation ( |l3| ) one can construct singular vector- 
potential A~ through the periodic solution A', A s = A' — 
V0' whereW 



{Vi ~ V) e x + {% - Xj)e y 
(x - Xi) 2 + {y- Vi) 2 



(18) 



Sum converges very slowly and the rate of convergence 
can be improved by Ewald construction 

V(ff = B ^ -Qv^ + Q^v exp ( _ 7 q 2) sin (Qr)+ 



E 



Q#0 

(yi - y)e x + (x - Xj)e y 
(x - x l ) 2 + {y- yi) 2 



■ exp 



(r - r,) 2 



4" 



(19) 



Parameter 7 is in principle arbitrary; we choose 7 = tt to 
assure good convergence for both sums. 

Self-consistency equations are solved by iteration. For 
a given vector-potential A' and gap function 'J, Eq. (||) 
is solved numerically. We start the iteration procedure 
with A' = and = const. Functions a and b are then 
used to calculate new values of A' and \& with help of self- 
consistency equations, and the whole process is repeated 
until the self-consistency is achieved. Equation (||) is 
also solved by iteration. We solved it in Fourier space 
where unknown are coefficients a. Iterative procedure 
for unknown coefficients a is 



o(Q) 



|*|(1 - a 2 )/2 - iuA s a] 



FT 



iQu 



(20) 



Maximum summation frequency was set as lo c = 2{)t 
(huj c — 207rfcsT c in physical units). Vortex lattice unit 
cell is divided into N x N mesh which turns into approx- 
imating auxiliary functions a and b by sum over N x N 
harmonics. We set N — 128. Integration over polar an- 
gle ip in each quadrant is performed by 30-points Gauss- 
Legendre integration. Iteration procedure fails for ex- 
tremely small lo, i.e. for very low temperatures t < 0.05. 



In all other cases it is recommended to damp the itera- 
tion process via a' new = {c x a new + c 2 a old )/(c 1 + c 2 ) with 
ci < C2- To manipulate with eqs. ([IJ) and, (|l5|), fast 
Fourier transform algorithm should be used.£3 

Now we concentrate on type-II superconductors with 
k ~ 1. For a particular value of k and t set of Eilcn- 
berger equations are solved numerically for various values 
of magnetic induction B, so we can construct free energy 
density F(B). However, in the experiment we actually 
control applied field H and not induction B. For a given 
applied filed H proper thermodynamic potential is Gibbs 
energy G(B,T) which is at minimum when the thermo- 
dynamic equilibrium is achieved, dG/dB = 0. This con- 
dition is rewritten via Hclmholtz free energy F(B,T) 



1 dF 
2P ~qb~- 



(21) 



If there are meta-stable states, as we expect in type-IIa 
superconductors, function B(H) is multi-valued and to 
find a true magnetic induction one should look for global 
minimum of Gibbs energy. In Fig. [I] we plotted typical 
example of magnetization curve with meta-stable states, 
i.e. for some values of H there are three different values 
of magnetic induction that correspond to the extrema of 
Gibbs energy. To find B which corresponds to the mini- 
mum energy one has to construct Gibbs energy density: 



G 



F - 2k 2 BH. 



(22) 



For example, we found that the minimum Gibbs energy 
for H = 0.34 correspond to the highest value of induction 
B = 0.20. Therefore from the behavior of Gibbs energy 
we can reconstruct the B(H) curve. From H — up to 
H = magnetic induction is zero and we have Meiss- 
ner state. Note that H*, is less than H c \ = lim H(B). 

At H = iJ*L there is first order phase transition from the 
Meissner state B = to the mixed state B = Bq (jump in 
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FIG. 3: Temperature dependence of magnetic field in the cen- 
ter of the vortex for parameters given in the text. 



magnetization). It means that energy of well separated 
vortices is bigger than vortex-lattice energy at B = B . 
due to attraction of vortices. For a fixed temperature 
the value of induction jump Bq decreases with increasing 
k and at some critical value k c , B q = 0. Then for all 
k > k c we have a second order phase transition from the 
Meissner to mixed state, i.e. type-lib superconductivity. 

Calculations are performed for various temperatures 
t = 0.2 — 0.9 to find the critical parameter k c that sep- 
arates type-IIa and type-lib superconductivity. The re- 
sults are shown on Fig. |[ Boundary line is qualitatively 
similar to the 3D case obtained by different numerical 
scheme. Boundary between type-I and type-II supercon- 
ductivity can be obtained from the condition H c = H C 2. 
Calculation itself is trivial and we did it merely to show 
the quantitative difference between 3D (shown by dashed 
line in Fig. ||) and 2D case. Compared to 3D case, bound- 
ary between type-I and type-II superconductivity in 2D 
superconductors is almost constant and changes only few 
percent down to T = 0. In 3D case there is possibility 
for cross-over from type-I to type-II superconductivity as 
one goes from T c to T = 0. 

We obtained selfconsistent solution for the order pa- 
rameter ^>(r,t) and field distribution A(r,t) for various 



k, B and temperature T. Only one interesting aspect 
will be presented: magnetic field profile as viewed via 
nuclear magnetic resonance (NMR) and muon spin rota- 
tion (jtiSR). These are standard experimental techniques 
that probe the magnetic field profile in the mixed state. 
NMR absorption intensity in a small frequency interval 
Aui of the ac field is proportional to the fraction of the 
VL cell area where the field value is in the correspond- 
ing interval Ato cx AH, i.e. signal is proportional to 
J 8(H — H(r))Ar. NMR line shape has the logarithmic 
singularity at the saddle point of magnetic field. In s- 
wave superconductor from the relative position of the 
NMR peak one can make conclusion about the VL struc- 
ture, j— . 

Kungtij measured NMR profile in V foil as a function of 
temperature. The applied field was perpendicular to the 
foil plane. In this geometry magnetic induction was fixed 
and almost temperature independent. It was found that 
the maximum magnetic field (field at the vortex center) 
is linear in T down to the lowest attainable temperature 
Q.1T C in that experiment. Delrieuo calculated magnetic 
field profile in the limit T = and small magnetiza- 
tion. He found singularities at the positions of field max- 
imum and minimum which have conical shape compared 
to parabolic shape near T c . It was shown that the slope 
dH max /dT is finite at T = 0. 

We calculated temperature dependence of H max in su- 
perconductor with k = 0.8 for fixed magnetic induction 
B = 0.16. The result is shown in Fig. ||. Overall be- 
havior H max (T) resemble the linear behavior reported 
inlI3. The result is in agreement with the calculation of 
Delrieu if the curve is extrapolated to T = 0. We note 
that in Ref.B H max has a more pronounced convex cur- 
vature. Therefore, in spite the value of magnetic pene- 
tration depth is saturated already below 0.4T C , magnetic 
field distribution is temperature dependent down to very 
low temperatures. 

In summary, we presented the new method to solve 
sclfconsistently Eilenberger equations in the mixed state. 
Its efficiency is demonstrated in the problem of 1ow-k 2D 
superconductors. Being faster then previous techniques, 
we hope it will stimulate further numerical studies of 
type-II superconductors. 
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